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Abstract 

In this file, we apply the nonequilibrium Green's function method to calculate the coefficient 
of thermal expansion (CTE) in multiple layers of graphene. We focus on the effect from different 
layer number N. 

Our main prediction is that: the increase of N can either enhance or weaken CTE, depending 
on the strength of the substrate interaction 7. If 7 < e, where e is the inter-layer interaction, the 
CTE will increase with increasing N. Otherwise, if 7 > e, CTE will decrease with increasing N. 

PACS numbers: 65.80.+n, 05.70.-a, 61.46.-w, 62.23.Kn 
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I. NONEQUILIBRIUM GREEN'S FUNCTION METHOD AND INTERACTION 
POTENTIAL 

We calculate the CTE of the multiple layers of graphene (MLG) using the nonequilibrium 
Green's function method, which includes all phonon modes automatically and can be applied 
to all temperatures as it is a pure quantum approach.- In this method, firstly, the averaged 
vibrational displacement of atom j, (uj), is calculated through: (uj) = ihGj, where Gj is 
the one point Green's function and can be calculated from its Feynman diagram expansion 
in terms of the nonlinear interaction. We consider the nonlinear interaction as, 

k, h 

H n = }^ -^~ U l U mUn + }^ -^—UoUpUqUr, (1) 
Imn opqr 4 

where u\ is the vibrational displacement of atom j, multiplied by the square root of its mass. 
The nonlinear force constant is extracted from "General Utility Lattice Program" (GULP) 2 
by the finite difference method. Then the CTE can be calculated explicitly by its definition, 

d(uj) 1 

a 3 = ^-f x -, 2 
dl Xj 

where Xj is the position of atom j along the expanding direction in the graphene sheet. We 
can get the final value of CTE by averaging ctj over all atoms. To avoid possible boundary 
effects, we have dropped 10% atoms which are on the two boundaries along the expanding 
direction. We apply periodic boundary condition in the perpendicular direction. Throughout 
this file, the z axis is perpendicular to the graphene sheet, and x y axes are lying in the 
graphene plane. 

The Tersoff potential is used to describe the intra-layer carbon carbon interactions. Its 
accuracy has been confirmed by comparison with results from ab initio density functional 
theory- For the inter-layer interaction, since the distance between two layers is out of the 
interaction range of Tersoff, we introduce Lennard- Jones potential, V(r) = 4e((a/r) 12 — 
(cr/r) 6 ), with e = 2.5 mev and a = 3.37 A. The value of the length parameter a is fitted 
to the space between adjacent layers in three-dimensional graphite 4 as 3.35 A. The energy 
parameter e is fitted to the phonon dispersion along TA direction in graphite as shown in 
Fig. [TJ which is comparable with the experimental results in Ref. |g. The cutoff for Lend- 
Jones potential is chosen as 20 A, which is large enough, i.e., all results do not change with 
further increasing this cutoff value. 
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FIG. 1: Phonon dispersion of three-dimensional graphite along TA direction in the Brillouin zone. 



The interaction between substrate and carbon atoms in the first graphene layer is simu- 
lated by the onsite potential in the z direction 7 : V = (7/2)1^ with 7 as the interaction force 
constant in unit of eV/ (A 2 u). u z is the vibrational displacement in the z direction multiplied 
by the square root of mass of the atom. There are three energy scales in the MLG system. 
The first one is the intra-layer interaction with the force constant in the order of ki n t ra ~ 1 
eV/ (A 2 u) which is estimated from Ref. |4|. In the following context we may simply say that 
the intra-layer interaction is in the order of 1.0. The second one is the inter-layer interaction 
with the force constant in the order of ki nter ~ 10~ 2 eV/(A 2 u) which is deduced from the 
above Lennard- Jones potential. Similarly, without mentioning we would say the inter-layer 
interaction e is in the order of 1CT 2 . The third one is the substrate interaction, which is 
described by the parameter 7. When we compare the strength of different interactions, we 
actually compare their corresponding force constants in the unit of eV/(A 2 u). 
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II. SIZE AND SUBSTRATE EFFECT 




T(K) 



FIG. 2: CTE v.s. T in single layer graphene sheet with different size, (length, width), and different 
substrate interaction 7. 

Now we investigate the CTE in single layer graphene sheet with different sizes and sub- 
strate interactions. Fig. [5] (a) shows that without substrate interaction, the CTE is very 
sensitive to the size of the system. The value of CTE decreases very quickly with the in- 
crease of length. For graphene with length larger than 40 A, the CTE is a very large negative 
number, -130 (10 _6 K _1 ), which means that the system with this size is thermally unstable. 
However, if the graphene sheet is put on a substrate, the situation is quite different. In case 
of weak substrate interaction, i.e. 7 < e as shown in figure (b), the size effect is greatly 
reduced. Yet it still have some important effect. If the substrate interaction is strong, i.e., 
7 > e as shown in figure (c), now the size effect is very small. If 7 3> e, CTE of all samples 
are very close to each other at all temperatures and positive in whole temperature range as 
shown in figure (d). 
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FIG. 3: Contribution to CTE from six lowest-frequency phonon modes in graphene sheet without 
substrate interaction. Blue solid line is for graphene with (length, width)=(10, 8.5) A, green dashed 
line is for (20, 17) A and red dotted line is for (20, 8.5) A. Insets are the vibrational morphology 
for the corresponding mode, (a), (b), (c) and (e) are the first four bending modes, (d) is a tearing 
mode, (f) is the longitudinal vibrational mode. 
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FIG. 4: Contribution to CTE from six lowest-frequency phonon modes in graphene sheet with 
substrate interaction 7 = 0.001. 

To understand the size and substrate effect on the CTE, we study the contribution to 
CTE from different phonon modes. In Fig. EJ we show contribution from the six lowest- 
frequency phonon modes. There is no substrate interaction in this figure, i.e., 7 = 0. The 
inset in each panel is the corresponding vibrational morphology of the system in this phonon 
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mode, (a), (b), (c) and (e) are the first four bending modes, (d) is an interesting tearing 
mode, which may be important in other thermal mechanical process. Among all six phonon 
modes, we find that the first bending mode shown in (a) dominants the value of CTE in 
all systems. Its contribution can reach as much as 90%. Due to the morphology of bending 
movement, it will induce contraction effect in the graphene sheet. So it leads to negative 
value of CTE. Since all of the first five phonon modes have negative effect on CTE while 
only the last mode has positive effect, the value of CTE is negative in all samples in case 
of no substrate interaction. Also, it is easier for longer system to bend,- so the negative 
effect from the bending mode increases rapidly with increasing length. That is the reason 
for more serious thermal contraction effect in longer system. 

If the substrate interaction is nonzero as shown in Fig. HI the value of the CTE is enhanced 
obviously and the difference between different sizes is narrowed. Now the contribution 
from the second bending mode is also very important. We can also see that the substrate 
interaction is more important in larger system, while it is less important in smaller system. 
Because the bending movement in larger graphene is more serious than the smaller system. 
When the substrate interaction is very strong, the graphene can not bend any more. So 
all bending modes do not contribute and only the sixth mode makes positive contribution 
to CTE. As a result, the CTE is positive in whole temperature range, and the difference 
between systems with different sizes is pretty small. 

III. N EFFECT 

Now we study the CTE in MLG with different layer numbers N (from 1 to 4) as shown in 
Fig. [5j In case of no substrate interaction (see panel (a)), CTE increases with the increase 
of N. Typically, when N changes from 1 to 2, the value of CTE shows a big jump. Because 
the inter-layer interaction makes the bending movement more difficult. So more graphene 
layers getting together can increase the thermal stability of the MLG if there is no substrate 
interaction. When the substrate interaction is not very strong, 7 < e, we can see similar 
but weaker effect, i.e., the CTE increases gradually with increasing N. If the substrate 
interaction is strong, 7 > e, we observe opposite phenomenon where the CTE decreases with 
increasing layer number N. At first glance, this is quite strange. However, it is actually 
physically understandable. If a piece of graphene is transferred onto a substrate with strong 
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FIG. 5: CTE in MLG with different layer number N = 1,2, 3, 4. The size for each layer is (length, 
width)=(10, 8.5) A. From (a) to (d) the strength of substrate interaction increases gradually. 

interaction, this graphene sheet is very stable. So the bending effect is small, leading to large 
value of CTE. When the layer number increases from 1 to 2, now the first graphene layer is 
very stable, but the second layer, which interacts with first layer with e, is not stable. As a 
result the thermal stability of the whole MLG decreases compared with N = 1. Similarly, 
the stability decreases gradually with further increasing N. When the substrate interaction 
is extremely large as shown in (d), we see similar and stronger effect. Fig. [6] shows the layer 
number dependence of CTE with more values of N. We can see that at all temperatures in 
the figure the CTE increases with increasing N in case of weak or no substrate interaction, 
and will decrease with increasing N in case of strong substrate interaction. The layer number 
dependence is more significant before N = 4. After N > 5, the difference is very small and 
the CTE reaches a saturate value independent of the substrate interaction. 

So our theoretical prediction is that: if the substrate interaction is weak or zero (j < e), 
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FIG. 6: CTE in MLG with layer number N from 1 to 10 at different temperatures. 

the value of CTE will increase with the increase of N; if the substrate interaction is strong 
(7 > e), the value of CTE will decrease with the increase of N. 

IV. POSSIBLE EXPERIMENTS 




L=10 




L=20 



FIG. 7: Clamp the graphene sample on the left boundary, and increase the size of the graphene. 
The lengths are 10, 20, 40 A in three figures. The large sample will be thermally unstable. 



10 





FIG. 8: Determine the relation between substrate interaction 7 and the inter-layer interaction e. 
(a) is a MLG sample with JV = 3 on a substrate. Apply a force F to the out-most graphene sheet 
to pull out the sample. If 7 < e, the sample will be peeled off the substrate as a whole completely 
as shown in (b). Otherwise, if 7 > e, the inner graphene sheet still stick on the substrate while the 
outer layers of the sample are peeled off as shown in (c). 
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In this section, we suppose two possible experiments to confirm our theoretical predic- 
tions. The first experiment is to examine the thermal instability of the graphene sheet 
as displayed in Fig. [7J Single layer graphene sheet with different sizes are picked out and 
clamped on the left boundary. It is understandable that the graphene sheet will be thermally 
unstable if the size is very large. According to our theoretical calculation, the typical length 
is about 40 A. 

The second experiment is to examine the layer number dependence of CTE. The MLG 
samples with different layer numbers N are transferred onto a substrate and clamped on 
the left boundary. Then the value of CTE in samples with different N can be measured by 
Lau's method.^ We mentioned that although the size of the system in our theoretical study 
is much smaller than the practical experimental samples, the physical mechanism should be 
the same in the experiment. According to our prediction, with the increase of N, the value 
of CTE will increase if 7 < e, or decrease if 7 > e. So it is very important to determine 
the relationship between 7 and e experimentally. As demonstrated by Fig. [SJ this can be 
done in a very simple way. Figure (a) is a MLG sample on the substrate, where a force F 
is applied onto the out-most graphene sheet, trying to pull out the MLG sample. If 7 < e, 
the whole MLG sample will be peeled off the substrate as shown in (b). On the other hand, 
if 7 > e, the first graphene layer will still stick to the substrate, while the other layers are 
peeled off the substrate shown in (c). 
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